#!/usr/bin/python
#
# this is alpha stage SW !
#
# DEPENDENCIES
# - numpy (install via package manager)
# - PILLOW (install via package manager)
#
#How to use it:
#python even-posterizer.py -Q img_6524.jpg
# all arguments not starting with "-" are considered images
#(intended for python2.7 and 3.x)
#
#any number of image files (wildcards allowed)
#(I presume you are working on linux)
# use "-h" for more info


from __future__ import division
from random import random
from sys import argv
#import math
import CS_convertor
import numpy as np
from os.path import split,splitext
from time import time


newsize=0  #0 for no resize
halfmode=False  # Leave part of image unchanged
withstats=False   #insert statistics on image
passresiduals=False  #Dithering
passressratio=0.5 #0-1, residuals are multiplied by this number
#forcebw=True  #force 0 for first value and 1 for last value (or 0 and 255 respectively)
blurradius=0
bandsize=40
bandmaxprc=25 # max % of sample pixels for single band
docleaning=False
squarecleaning=False
bandmultipl=1.0
sqinnersize=2
outlining=False  #create image with outlines
mpart=""
sqpart=""

version="0.1.7"

def insert_text(x,y,text,draw):
	draw.text((x+1,y-1),text,(0,0,0))	
	draw.text((x+1,y+1),text,(0,0,0))	
	draw.text((x-1,y-1),text,(0,0,0))	
	draw.text((x-1,y+1),text,(0,0,0))
	draw.text((x ,y ),text,(255,255,0))

def get_minmax(negsamples,possamples):
	if len(negsamples)>0:
		minv=-negsamples[-1]
	else:
		minv=possamples[0]
	if len(possamples)>0:
		maxv=possamples[-1]
	else:
		maxv=-negsamples[0]
	return minv,maxv



def get_position(Value,Limitsarray):
	for i in range(len(Limitsarray)):
		if Value<=Limitsarray[i]:
			return i
	else:
		print("out of range!!!")
		print (Value)
		print (Limitsarray)
		exit(3)
		#return postlevel-1
		
def count_val_in_range(array,minvalue,maxvalue):
	count=0
	for item in (array):
		if item >= minvalue and item <maxvalue:
			count=count+1
	return count

def get_final(limitsarray):
	finalarray=[]
	for i in range(len(limitsarray)):
		if i==0:
			finalarray.append((0+limitsarray[i])/2.0)
		else:
			finalarray.append((limitsarray[i-1]+limitsarray[i])/2.0)	
	return finalarray


def get_limits(samplesarray,defaultstep,bandmaxprclocal,verbose=False):
	#splitting L axis into band and assigning final value for every band
	#limits contain upper value (last one is 100)
	startingvalue=0
	step=defaultstep
	finallimits=[]
	while True:
		if verbose:
			print ("starting value: "+str(startingvalue)+", step: "+str(step))
		count=count_val_in_range(samplesarray,startingvalue,startingvalue+step)
		if verbose:
			print (" Count? "+str(count))
		if step<0:
			break
		if verbose:
			print (count,samplescount*bandmaxprc/100.0)
		if count<samplescount*bandmaxprclocal/100.0 or step<=2:
			if startingvalue+step<100:
				finallimits.append(startingvalue+step)
				startingvalue=startingvalue+step
				step=bandsize			
			else:
				finallimits.append(100)
				break			

		else:
			if verbose:
				print ("reducing step")
			step=step-1
	return finallimits		


			
try:
	from PIL import Image      # necessary for manipulation with images
	from PIL import ImageFilter
	from PIL import ImageDraw,ImageFont
except:
	print (" Warning: Failed to load Pillow (Imaging) modules!")
	print (" Pillow (fork of PIL) should be available in repositories")
	print (" of your distribution, please install it first.")
	exit()


def cleaning_sq(innersize,workingarray):
	changed=0
	for x in range(0,imgX-innersize-2):
		for y in range(0,imgY-innersize-2):	
			minvalue=100
			maxvalue=-100
			#scanning outer pixels of square
			for i in range(innersize+1):
				#A
				startx=x
				starty=y
				#step=x+i
				localvalue=workingarray[startx+i,starty]
				if localvalue<minvalue:
					minvalue=localvalue
				if localvalue>maxvalue:
					maxvalue=localvalue
				#B
				startx=x
				starty=y
				#step=y+i
				localvalue=workingarray[startx,starty+i]
				if localvalue<minvalue:
					minvalue=localvalue
				if localvalue>maxvalue:
					maxvalue=localvalue
				#C				
				startx=x+innersize+1
				starty=y
				#step=y+i					
				localvalue=workingarray[startx,starty+i]
				if localvalue<minvalue:
					minvalue=localvalue
				if localvalue>maxvalue:
					maxvalue=localvalue			
				#D
				startx=x
				starty=y+innersize+1
				#step=x+i					
				localvalue=workingarray[startx+i,starty]
				if localvalue<minvalue:
					minvalue=localvalue
				if localvalue>maxvalue:
					maxvalue=localvalue		



			
			#changing inner pixels
			for x_l in range (1,1+innersize):
				for y_l in range (1,1+innersize):
					if workingarray[x+x_l,y+y_l]<minvalue:
						workingarray[x+x_l,y+y_l]=minvalue
						changed=changed+1
					if workingarray[x+x_l,y+y_l]>maxvalue:
						workingarray[x+x_l,y+y_l]=maxvalue
						changed=changed+1		
	return changed			
					
					
def cleaning2(othertresh,sameothertresh,valuediff,workingarray,mask):
	changed=0
	skipped=0
	for x in range(1,imgX-1):
		for y in range(1,imgY-1):
			if mask[x,y]==1:
				skipped=skipped+1
				continue
			nears={}
			newvalue=-1000
			othergoodcount=0
			othercount=0
			#othersamecount=0
			for relpos in nearpoints:
				nearvalue=workingarray[x+relpos[0],y+relpos[1]]
				if nearvalue==workingarray[x,y]:
					pass
				else:
					othercount=othercount+1
					if abs(nearvalue-workingarray[x,y])<valuediff:
						othergoodcount=othergoodcount+1
						if nearvalue in nears: #.has_key(nearvalueL):
							nears[nearvalue]=nears[nearvalue]+1
						else:
							nears[nearvalue]=1

			#statistics done
			if othercount==0:
				mask[x,y]=1 
				continue
			if othercount>=othertresh:
				#finding neighbourd with count same or above sameothertresh:
				goodcount=0
				newvalue=-1000
				for k,v in nears.items() :
					if v>goodcount:
						newvalue=k
						goodcount=v
				if goodcount>=sameothertresh:
					changed=changed+1
					#print (workingarray[x,y],newvalue)
					workingarray[x,y]=newvalue

					#changing status of near pixels
					for relpos in nearpoints:
						mask[x+relpos[0],y+relpos[1]]=0
	#print ("   skipped:" +str(skipped))
	return changed


def do_outlines(final):
	for x in range(1,imgX-1):
		for y in range(1,imgY-1):
			diff=0
			centralvalueL=Loutput[x,y]
			centralvalueA=Aoutput[x,y]
			centralvalueB=Boutput[x,y]
			for relpos in [[0,1], [1,1], [1,0],[1,-1],[0,-1],[-1,-1],[-1,0], [-1,1]]:
				nearvalueL=Loutput[x+relpos[0],y+relpos[1]]
				nearvalueA=Aoutput[x+relpos[0],y+relpos[1]]
				nearvalueB=Boutput[x+relpos[0],y+relpos[1]]
				tmpdiff=5.0 * (abs(nearvalueL-centralvalueL) + abs(centralvalueA-nearvalueA) + abs(centralvalueB-nearvalueB))
				if diff<tmpdiff:
					diff=tmpdiff
				#diff=diff+abs(nearvalueL-centralvalueL) + abs(centralvalueA-nearvalueA) + abs(centralvalueB-nearvalueB)
			#print diff
			final[x,y]=diff
				
				

def cleaning(treshold,difftresh):
	lchanged=0
	achanged=0
	bchanged=0
	for x in range(1,imgX-1):
		for y in range(1,imgY-1):	
			nearsL={}
			newvalueL=-1000
			nearsA={}
			newvalueA=-1000
			nearsB={}
			newvalueB=-1000
			undoablecountL=0
			undoablecountA=0
			undoablecountB=0
			for relpos in [[0,1], [1,1], [1,0],[1,-1],[0,-1],[-1,-1],[-1,0], [-1,1]]:
				nearvalueL=Loutput[x+relpos[0],y+relpos[1]]
				nearvalueA=Aoutput[x+relpos[0],y+relpos[1]]
				nearvalueB=Boutput[x+relpos[0],y+relpos[1]]
				if nearvalueL==Loutput[x,y]:
					undoablecountL=undoablecountL+1
				elif abs(nearvalueL-Loutput[x,y])>difftresh:
					undoablecountL=undoablecountL+1					
				elif nearvalueL in nearsL: #.has_key(nearvalueL):
					nearsL[nearvalueL]=nearsL[nearvalueL]+1
				else:
					nearsL[nearvalueL]=1
				if nearvalueA==Aoutput[x,y]:
					undoablecountA=undoablecountA+1
				elif abs(nearvalueA-Aoutput[x,y])>difftresh:
					undoablecountA=undoablecountA+1	
				elif nearvalueA in nearsA: #.has_key(nearvalueA):
					nearsA[nearvalueA]=nearsA[nearvalueA]+1
				else:
					nearsA[nearvalueA]=1
				if nearvalueB==Boutput[x,y]:
					undoablecountB=undoablecountB+1
				elif abs(nearvalueB-Boutput[x,y])>difftresh:
					undoablecountB=undoablecountB+1	
				elif nearvalueB in nearsB: #.has_key(nearvalueB):
					nearsB[nearvalueB]=nearsB[nearvalueB]+1
				else:
					nearsB[nearvalueB]=1


			if undoablecountL<=8-treshold:
				goodcount=0
				for k,v in nearsL.items() :
					if v>goodcount:
						newvalueL=k
						goodcount=v
			if undoablecountA<=8-treshold:
				goodcount=0
				for k,v in nearsA.items() :
					if v>goodcount:
						newvalueA=k
						goodcount=v
			if undoablecountA<=8-treshold:
				goodcount=0
				for k,v in nearsB.items() :
					if v>goodcount:
						newvalueB=k
						goodcount=v


			if newvalueL>-1000:
				lchanged=lchanged+1
				Loutput[x,y]=newvalueL
				#print (newvalue)
			if newvalueA>-1000:
				achanged=achanged+1
				Aoutput[x,y]=newvalueA
			if newvalueB>-1000:
				bchanged=bchanged+1
				Boutput[x,y]=newvalueB

	print (lchanged,achanged,bchanged)

print ("[even-cartoonizer.py, version: "+str(version)+", use '-h' for help and bit more info]")

files=[]
#parsing arguments
for arg in argv[1:]:
	if arg.startswith("-b"):
		#print (arg[2:])
		blurradius=1
		print ("Bluring ...")
	if arg=="-o":	
		outlining=True
		print (" Creating also img with outlines...")
	if arg=="-H":	
		halfmode=True
		print (" Halfmode...")
	if arg.startswith("-Q"):
		squarecleaning=True
		if arg=="-Q":	
			print (" Cleaning (square mode)...")
			sqpart="_Q"
		#print (arg[2:])
		else:
			try:
				sqinnersize=int(arg[2:])
			except:
				print ("Failed to get inner size for sq cleaning:" + str(arg[2:])+" ???")
				sqinnersize=0
			
			if (sqinnersize<1):
				sqinnersize=1
			if (sqinnersize>20):
				sqinnersize=20
			print (" New ineer size for square cleaning: "+str(sqinnersize))
			sqpart="_Q"+str(sqinnersize)
	if arg.startswith("-m"):

		try:
			temp=float(arg[2:])
		except:
			print (" ! Failed to get banding multiplicator:" + str(arg[2:])+" ??? (0.01-3)")
			continue
			
		if (temp<0.01 or temp>3.0):
			print (" Banding multiplicator out of range (0.01-3) :" + str(temp))
			continue
		print (" Banding multiplicator: "+str(temp))
		bandsize=bandsize*temp
		bandmaxprc=bandmaxprc*temp
		mpart="_m"+str(temp)[:3]
		bandmultipl=temp

	if arg=="-C":	
		docleaning=True
		print (" Cleaning on...")
	if arg=="-S":	
		withstats=True
		print (" Statistics will be printed on image...")
	if arg.startswith("-s"):
		#print (arg[2:])
		try:
			newsize=int(arg[2:])
		except:
			print ("Failed to get new size:" + str(arg[2:])+" ???")
			newsize=0
		
		if (newsize<10):
			newsize=10
		if (newsize>2500):
			newsize=2500
		print (" New size: "+str(newsize))
	if arg.startswith("-h") or arg.startswith("--help"):
		print (" ")
		print ("even-cartoonizer.py "+version)
		print ("Available switches")
		print ("-s#  where # is width of image (no space)")	
		print ("-b  blur the image before processing(no intensity/number is to be given)")	
		print ("-Q#  Square-cleaning mode with inner size of the square (without number the default is used)")	
		print ("-m#  Band width multiplicator. Range: 0.01-3. If >1 colors count is reduced, if <1 colors count is increased")	
		print ("-H  halfmode")	
		print ("-o  create also separate image with lines")	
		print ("-S	put statistics on the image")
		print ("All switches not starting with '-' are considered images")
		print (" ")
		print ("This is prototype/alpha level software... ")	
		print ("Contact: tiborb95 at gmail dot com ")		
		exit (0)
		
	if not arg.startswith("-"):
		files.append(arg)



if len(files)==0:
	print ("No files/images given...")

gammedValues=[ (value/255.0) ** (1/2.2) for value in range(256) ]	

imgcount=0	
for image in files:
	imgcount=imgcount+1	
	
	# trying to open image
	try:
		sourceImg = Image.open(image)
		sourceImg_loaded=sourceImg.load()
	except:
		print ("Cannot open "+image+" image, skipping this one...")
		continue
	
	print ("\nProcessing: "+image+"      "+str(imgcount)+"/"+str(len(files)))
	
	corename=splitext(split(image)[1])[0]	
	
	#resizing if needed
	imgX,imgY=sourceImg.size
	if newsize>0 and imgX > newsize:
		print (" Resizing...")

		sourceImg=sourceImg.resize((newsize,int(newsize/imgX*imgY)),Image.ANTIALIAS)
		if blurradius>0:
			print (" Blurring...")
			sourceImg = sourceImg.filter(ImageFilter.GaussianBlur)
		sourceImg_loaded=sourceImg.load()
	
	
	# BODY
	imgX,imgY=sourceImg.size
	
	Lsamples=[]
	Anegsamples=[]
	Bnegsamples=[]
	Apossamples=[]
	Bpossamples=[]
	samplescount=0
	
	for x in range(0,imgX,int(imgX/20)):
		for y in range(0,imgY,int(imgY/20)):
			R=sourceImg_loaded[x,y][0]
			G=sourceImg_loaded[x,y][1]
			B=sourceImg_loaded[x,y][2]
			l,a,b = CS_convertor.rgb_to_cielab(R,G,B)
			Lsamples.append(l)
			if a<-0.5:
				Anegsamples.append(-a)
			if a>0.5:
				Apossamples.append(a)			
			if b<-0.5:
				Bnegsamples.append(-b)
			if b>0.5:
				Bpossamples.append(b)	
			samplescount=samplescount+1

		
			
	# ordering the samples lists
	Lsamples.sort()
	Anegsamples.sort()
	Bnegsamples.sort()
	Apossamples.sort()
	Bpossamples.sort()
	#print (Anegsamples)
	#print (len(Anegsamples),len(Bnegsamples), len(Apossamples), len(Bpossamples))
	#printing samples range
	print (" CH1 range: "+str(Lsamples[0])[:6]+"  -  "+str(Lsamples[-1])[:6])
	minv,maxv=get_minmax(Anegsamples,Apossamples)
	print (" CH2 range: "+str(minv)[:6]+"  -  "+str(maxv)[:6])
	minv,maxv=get_minmax(Bnegsamples,Bpossamples)
	print (" CH3 range: "+str(minv)[:6]+"  -  "+str(maxv)[:6])
	

	
	Llimits=get_limits(Lsamples,bandsize,bandmaxprc,False)
	Lfinal=get_final(Llimits)
	Aneglimits=get_limits(Anegsamples,bandsize,bandmaxprc,False)
	Anegfinal=get_final(Aneglimits)
	Aposlimits=get_limits(Apossamples,bandsize,bandmaxprc,False)
	Aposfinal=get_final(Aposlimits)
	Bneglimits=get_limits(Bnegsamples,bandsize,bandmaxprc,False)
	Bnegfinal=get_final(Bneglimits)
	Bposlimits=get_limits(Bpossamples,bandsize,bandmaxprc,False)
	Bposfinal=get_final(Bposlimits)
	print (" Count of final colors per channel: L: "+str(len(Lfinal))+", A: "+str(len(Anegfinal)+len(Aposfinal)+1)+", B: "+str(len(Bnegfinal)+len(Bposfinal)+1))


	# creating new image
	newImage = Image.new( 'RGB', (imgX,imgY), "black") #
	#creating empty arrays for LAB 
	Loutput=np.zeros((imgX,imgY))
	Aoutput=np.zeros((imgX,imgY))
	Boutput=np.zeros((imgX,imgY))
	
	for x in range(0,imgX):
		if ((x+1)%200==0):
			print (" column: "+str(x+1)+"/"+str(imgX))
		rres=0
		gres=0
		bres=0
		for y in range(0,imgY):
			R=sourceImg_loaded[x,y][0]
			G=sourceImg_loaded[x,y][1]
			B=sourceImg_loaded[x,y][2]
			
			#if halfmode and x<imgX/2.5:
				#newImage.putpixel((x,y),(R,G,B))
				#continue
			#elif halfmode and x<imgX/2.5+5:
				#newImage.putpixel((x,y),(255-R,255-G,255-B))
				#continue

			l,a,b=CS_convertor.rgb_to_cielab(R,G,B)
			
			if passresiduals:
				lres=rres*passressratio
				ares=gres*passressratio
				bres=bres*passressratio
				if ((x+y*2)%7==0):
					#print (" no res passing")
					pass
				else:
					#print (" res passing")
					l=l-lres
					a=a-ares
					b=b-bres			

			lnew=Lfinal[get_position(l,Llimits)]
			#print (l,lnew,get_position(l,Llimits))
			if a>0.5:
				anew=Aposfinal[get_position(a,Aposlimits)]
			elif a<-0.5:
				anew=-Anegfinal[get_position(-a,Aneglimits)]
			else:
				anew=0
			if b>0.5:
				bnew=Bposfinal[get_position(b,Bposlimits)]
			elif b<-0.5:
				bnew=-Bnegfinal[get_position(-b,Bneglimits)]
			else:
				bnew=0
			#print (b,bnew)


					
			if passresiduals:
				lres=lnew-l					
				ares=anew-a	
				bres=bnew-b	


			Loutput[x,y]=lnew
			Aoutput[x,y]=anew
			Boutput[x,y]=bnew

	nearpoints=[[0,1], [1,1], [1,0],[1,-1],[0,-1],[-1,-1],[-1,0], [-1,1]]
	#square cleaning
	if squarecleaning:
		print (" Changed pixels per channel   :      L      A      B   Overall            sec")
		cleaned_sum=1000000
		while cleaned_sum>imgX*imgY*0.005:
			t1 = time()
			lchanged=cleaning_sq(sqinnersize,Loutput) #other values, same values in other, difftreshols
			achanged=cleaning_sq(sqinnersize,Aoutput) 
			bchanged=cleaning_sq(sqinnersize,Boutput) 
			t2 = time()
			cleaned=lchanged+achanged+bchanged
			print (" squarecleaning               : "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(round(t2-t1,1))[:5].rjust(6))
		
			t1 = time()
			Lmask=np.zeros((imgX,imgY))
			Amask=np.zeros((imgX,imgY))
			Bmask=np.zeros((imgX,imgY))
			lchanged=cleaning2(5,5,35,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(5,5,35,Aoutput,Amask) 
			bchanged=cleaning2(5,5,35,Boutput,Bmask) 
			t2 = time()
			cleaned1=lchanged+achanged+bchanged
			print (" standard cleaning (pass 1)   : "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned1/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(round(t2-t1,1))[:5].rjust(6))

			t1 = time()
			lchanged=cleaning2(5,5,35,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(5,5,35,Aoutput,Amask) 
			bchanged=cleaning2(5,5,35,Boutput,Bmask) 
			t2 = time()
			cleaned2=lchanged+achanged+bchanged
			print (" standard cleaning (pass 2)   : "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned2/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(round(t2-t1,1))[:5].rjust(6))
			
				
			cleaned_sum=cleaned+cleaned1+cleaned2



	#cleaning
	Lmask=np.zeros((imgX,imgY))
	Amask=np.zeros((imgX,imgY))
	Bmask=np.zeros((imgX,imgY))

	if docleaning:
		
		#single pixel cleaning:
		t1 = time()
		lchanged=cleaning2(8,0,100,Loutput,Lmask) #other values, same values in other, difftreshols
		achanged=cleaning2(8,0,100,Aoutput,Amask) 
		bchanged=cleaning2(8,0,100,Boutput,Bmask) 
		t2 = time()
		cleaned1=lchanged+achanged+bchanged
		print (" single px cleaned up: "+str(lchanged)+" "+str(achanged)+" "+str(bchanged)+" = "+str(cleaned1/imgX/imgY*100.0)[:5]+"% , time: "+str(t2-t1))
			
		
		cleaned1=1000000
		while cleaned1>imgX*imgY*0.002:
			t1 = time()
			lchanged=cleaning2(5,5,25,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(5,5,25,Aoutput,Amask) 
			bchanged=cleaning2(5,5,25,Boutput,Bmask) 
			t2 = time()
			cleaned1=lchanged+achanged+bchanged
			print (" cleaned up pixels (pass 2): "+str(lchanged)+" "+str(achanged)+" "+str(bchanged)+" = "+str(cleaned1/imgX/imgY*100.0)[:5]+"% , time: "+str(t2-t1))
			
		
		cleaned1=1000000
		while cleaned1>imgX*imgY*0.002:
			t1 = time()
			lchanged=cleaning2(6,5,50,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(6,5,50,Aoutput,Amask) 
			bchanged=cleaning2(6,5,50,Boutput,Bmask) 
			t2 = time()
			cleaned1=lchanged+achanged+bchanged
			print (" cleaned up pixels (pass 3): "+str(lchanged)+" "+str(achanged)+" "+str(bchanged)+" = "+str(cleaned1/imgX/imgY*100.0)[:5]+"% , time: "+str(t2-t1))

		#single pixel cleaning:
		t1 = time()
		lchanged=cleaning2(8,0,100,Loutput,Lmask) #other values, same values in other, difftreshols
		achanged=cleaning2(8,0,100,Aoutput,Amask) 
		bchanged=cleaning2(8,0,100,Boutput,Bmask) 
		t2 = time()
		cleaned1=lchanged+achanged+bchanged
		print (" single px cleaned up: "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned1/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(t2-t1)[:5].rjust(6))
		
		cleaned1=1000000
		while cleaned1>imgX*imgY*0.0005:
			t1 = time()
			lchanged=cleaning2(5,5,25,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(5,5,25,Aoutput,Amask) 
			bchanged=cleaning2(5,5,25,Boutput,Bmask) 
			t2 = time()
			cleaned1=lchanged+achanged+bchanged
			print (" cleaned up pixels (pass 5): "+str(lchanged)+" "+str(achanged)+" "+str(bchanged)+" = "+str(cleaned1/imgX/imgY*100.0)[:5]+"% , time: "+str(t2-t1))
			
		
		cleaned1=1000000
		while cleaned1>imgX*imgY*0.0005:
			t1 = time()
			lchanged=cleaning2(6,6,100,Loutput,Lmask) #other values, same values in other, difftreshols
			achanged=cleaning2(6,6,100,Aoutput,Amask) 
			bchanged=cleaning2(6,6,100,Boutput,Bmask) 
			t2 = time()
			cleaned1=lchanged+achanged+bchanged
			print (" cleaned up pixels (pass 6): "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned1/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(t2-t1)[:5].rjust(6))
	
		#single pixel cleaning:
		t1 = time()
		lchanged=cleaning2(8,0,100,Loutput,Lmask) #other values, same values in other, difftreshols
		achanged=cleaning2(8,0,100,Aoutput,Amask) 
		bchanged=cleaning2(8,0,100,Boutput,Bmask) 
		t2 = time()
		cleaned1=lchanged+achanged+bchanged
		print (" single px cleaned up: "+str(lchanged).rjust(6)+" "+str(achanged).rjust(6)+" "+str(bchanged).rjust(6)+" = "+str(cleaned1/imgX/imgY*100.0)[:5].rjust(6)+"% , time: "+str(t2-t1)[:5].rjust(6))
	
	
	if (outlining):
		Lines_array=np.zeros((imgX,imgY))
		do_outlines(Lines_array)
		linesImage = Image.new( 'RGB', (imgX,imgY), "black")
		for x in range(0,imgX):
			for y in range(0,imgY):	
				linesImage.putpixel((x,y),( int(255-Lines_array[x,y]) , int(255-Lines_array[x,y]) , int(255-Lines_array[x,y])))
			
	
	print (" Converting to RGB...")
	LAB_TO_RGB={}
	for x in range(0,imgX):
		for y in range(0,imgY):	
		
			if halfmode and x<imgX/2.5:
				newImage.putpixel((x,y),sourceImg_loaded[x,y])
				continue
			elif halfmode and x<imgX/2.5+1:
				newImage.putpixel((x,y),(255,255,255))
				continue
			elif halfmode and x<imgX/2.5+2:
				newImage.putpixel((x,y),(0,0,0))
				continue		
			hashableint=int(10349*Loutput[x,y] + 107*Aoutput[x,y] +Boutput[x,y] )
			if 	hashableint in LAB_TO_RGB: #.has_key(hashableint):
				Rnew,Gnew,Bnew=LAB_TO_RGB[hashableint]
			else:
				Rnew,Gnew,Bnew=CS_convertor.cielab_to_rgb(Loutput[x,y],Aoutput[x,y],Boutput[x,y])
				LAB_TO_RGB[hashableint]=[Rnew,Gnew,Bnew]
				#print ("  Calculating new Lab values for color with hash: "+str(hashableint)+" (doing: "+str(x)+"/"+str(imgX)+")")
				#print (Loutput[x,y],Aoutput[x,y],Boutput[x,y])
			#Rnew,Gnew,Bnew=CS_convertor.cielab_to_rgb(Loutput[x,y],Aoutput[x,y],Boutput[x,y])
			newImage.putpixel((x,y),( int(round(Rnew)),int(round(Gnew)),int(round(Bnew)) ) )
	
	#inserting text
	startypos=0.97
	if withstats:
		startypos=0.97
		draw = ImageDraw.Draw(newImage)
		if (squarecleaning):
			text="Square kernel cleaning on, ammount (inner size in px): "+str(sqinnersize)
			insert_text(int(imgX*0.1),int(imgY*startypos),text,draw)
			startypos=startypos-0.02
		if (passresiduals):
			text="Dithering on, weight="+str(passressratio)
			insert_text(int(imgX*0.1),int(imgY*startypos),text,draw)
			startypos=startypos-0.02
		if not (bandmultipl==1.0):
			text="Banding width multiplicator="+str(bandmultipl)
			insert_text(int(imgX*0.1),int(imgY*startypos),text,draw)
			startypos=startypos-0.02
		textR="L values: "+', '.join([str(numb)[:4] for numb in Lfinal])
		insert_text(int(imgX*0.1),int(imgY*(startypos-0.04)),textR,draw)
		textG="A values: "+', '.join([str(-numb)[:5] for numb in sorted(Anegfinal,reverse=True)]) + ", 0 ,"+ ', '.join([str(numb)[:4] for numb in Aposfinal])
		insert_text(int(imgX*0.1),int(imgY*(startypos-0.02)),textG,draw)
		textB="B values: "+', '.join([str(-numb)[:5] for numb in sorted(Bnegfinal,reverse=True)]) + ", 0 ,"+ ', '.join([str(numb)[:4] for numb in Bposfinal])
		insert_text(int(imgX*0.1),int(imgY*(startypos-0.00)),textB,draw)


	
		
	#saving new image
	if (passresiduals):
		dithpart="_D"+str(passressratio)
	else:
		dithpart=""
	newname="post_"+corename+dithpart+sqpart+mpart+".jpg"
	print( " Saving as: "+newname)
	newImage.save(newname,quality=95)
	newname="post_"+corename+dithpart+sqpart+mpart+".png"
	print( " Saving as: "+newname)
	newImage.save(newname,optimize=True)
	
	if (outlining):
		newname="post_"+corename+dithpart+sqpart+mpart+"_outlines.png"	
		print( " Saving lines as: "+newname)			
		linesImage.save(newname,optimize=True)
